{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this problem we tackle the learning problem for the canonical case of one-dimensional functions. This seemingly simple problem provides a rich test bed for understanding several interesting concepts both theoretically and empirically since almost everything can be visualized in this domain.\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "from ipywidgets import interactive\n",
    "import ipywidgets as widgets\n",
    "from ipywidgets import fixed\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part a): Generating training data\n",
    "\n",
    "Run the  cells in part a) (i.e upto the beginning of part b)) and visualize the noisy training data. Play around with the standard deviation of the additive white gaussian noise. What do you observe?\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "### Generating $x_i$\n",
    "Before we can learn a function from its samples we must first generate the training data $(x_i, y_i)_{i=1,\\dots,n}$. Here $n$ is the number of training samples. In this problem we consider two ways of sampling $x_i$ for training data. <br>\n",
    "1. $x_i$ sampled uniform randomly from $[-1,1]$.\n",
    "2. $x_i$ from an evenly spaced grid on the interval $[-1,1]$.\n",
    "   For example, for $n$ = 4, then we have the samples (-1, -0.5, 0, 0.5). Note that the endpoint 1, is not included in our training set.\n",
    "   This kind of evenly spaced samples gives rise to interesting properties of the feature matrix when using Fourier features and we will see this later in this problem.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def generate_x(n, x_type, x_low=-1, x_high=1):\n",
    "    if x_type == 'grid':\n",
    "        x = np.linspace(x_low, x_high, n, endpoint = False).astype(np.float64)\n",
    "\n",
    "    elif x_type == 'uniform_random':\n",
    "        x = np.sort(np.random.uniform(x_low, x_high, n).astype(np.float64))\n",
    "        #Note that for making it easy for plotting we sort the randomly sampled x in ascending order\n",
    "    else:\n",
    "        raise ValueError\n",
    "\n",
    "\n",
    "    return x\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Generating $y_i$\n",
    "Now that we have $x_i$  we can use it to generate $y_i = f(x_i)$ where $f$ is the true function that we are trying to learn.  In this problem we consider the following:\n",
    "1. $y = x$.\n",
    "2. $y = cos(2\\pi x)$.\n",
    "3. $y = sgn(x) =  \\begin{cases}\n",
    "         -1, & x \\leq 0 \\\\\n",
    "         1, & x = 0\n",
    "        \\end{cases}$\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def mysign(x):\n",
    "    y = np.sign(x)\n",
    "    y[x == 0] = 1\n",
    "    return y\n",
    "\n",
    "def generate_y(x, f_type):\n",
    "    if f_type == 'x':\n",
    "        y = x\n",
    "\n",
    "    elif f_type == 'cos2':\n",
    "        y = np.cos(2*np.pi * x)\n",
    "\n",
    "    elif f_type == 'sign':\n",
    "        y = mysign(x)\n",
    "\n",
    "    else:\n",
    "        raise ValueError\n",
    "\n",
    "    return y\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualizing training data\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_training_data(x_type):\n",
    "    n = 64\n",
    "    x_true = generate_x(x_type = 'grid', n=1000)\n",
    "    x_train = generate_x(x_type=x_type, n=n)\n",
    "    labels = ['y=x', 'y=cos(2${\\pi}$x)', 'y=sgn(x)']\n",
    "    for k,f_type in enumerate(['x', 'cos2', 'sign']):\n",
    "        y_true = generate_y(x=x_true, f_type=f_type)\n",
    "        y_train = generate_y(x=x_train, f_type=f_type)\n",
    "        plt.plot(x_true, y_true, linewidth = 0.5)\n",
    "        plt.ylabel('y')\n",
    "        plt.xlabel('x')\n",
    "        plt.scatter(x_train, y_train, marker='o', label = labels[k])\n",
    "    plt.legend(bbox_to_anchor  = (1.03, 0.97))\n",
    "    plt.show()\n",
    "\n",
    "\n",
    "\n",
    "slider = widgets.RadioButtons(\n",
    "    options=['uniform_random', 'grid'],\n",
    "    description='x_type:',\n",
    "    disabled=False\n",
    ")\n",
    "\n",
    "\n",
    "interactive_plot = interactive(plot_training_data, x_type=slider)\n",
    "output = interactive_plot.children[-1]\n",
    "interactive_plot\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Noise in training data\n",
    "\n",
    "Often in practice our training data is noisy. We model this in our problem by assuming that our samples for y are corrupted by Additive White Gaussian Noise (AWGN).  The true function is given by $y = f(x)$. The amount of noise is controlled by the standard deviation of the noise which we denote by awgn_std. The noiseless case corresponds to awgn_std = 0\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def add_awgn_noise(y, awgn_std=0):\n",
    "    noise = np.random.normal(0, awgn_std, y.shape)\n",
    "    y_noisy = y + noise\n",
    "    return y_noisy\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualizing noise in training data\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_noisy_training_data(awgn_std):\n",
    "    n = 64\n",
    "    np.random.seed(7)\n",
    "    x_true = generate_x(x_type = 'grid', n=1000)\n",
    "    x_train = generate_x(x_type='uniform_random', n=n)\n",
    "    f_type = 'x'\n",
    "    y_true = generate_y(x=x_true, f_type=f_type)\n",
    "    y_train_clean = generate_y(x=x_train, f_type=f_type)\n",
    "    y_train = add_awgn_noise(y_train_clean, awgn_std=awgn_std)\n",
    "    plt.plot(x_true, y_true, linewidth = 0.5, label = 'True function')\n",
    "    plt.ylabel('y')\n",
    "    plt.xlabel('x')\n",
    "    plt.ylim([-4,4])\n",
    "    plt.scatter(x_train, y_train, marker='o', label = 'Training samples')\n",
    "    plt.legend(loc = 'upper right', bbox_to_anchor  = (1.43, 0.97))\n",
    "\n",
    "    plt.show()\n",
    "\n",
    "\n",
    "slider= widgets.FloatLogSlider(\n",
    "    value=-3,\n",
    "    base=10,\n",
    "    min=-3, # max exponent of base\n",
    "    max=1, # min exponent of base\n",
    "    step=0.2, # exponent step\n",
    "    description='awgn_std',\n",
    "    continuous_update= False\n",
    ")\n",
    "\n",
    "interactive_plot = interactive(plot_noisy_training_data, awgn_std=slider)\n",
    "output = interactive_plot.children[-1]\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part b): Featurization- Lifting the training data\n",
    "\n",
    "\n",
    "Run the cell in part b) and understand visually how the polynomial features and Fourier features look like. Intuitively which choice of feature type will be easier to learn the function $y=x$? What about the function $y = cos (2 \\pi x)$?\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "We would like to learn the function from training samples by performing linear regression. However if the true function is not affine in $x$, we cannot hope to learn the function well. The key is to first lift the training data, i.e featurize it and express it in a richer space which we can then use to perform linear regression. For instance we could think of a polynomial featurization of $x$ given by $[1, x, x^2]$. Using this featurization if the true function were quadratic (say $y = x^2$),  by performing linear regression in the feature space we will succeed in learning a good approximation to the function. We will consider two types of featurizations that lift the one dimensional data into a d-dimensional feature space.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Polynomial features\n",
    "We consider the d-dimensional features given by the Vandermonde polynomials:\n",
    "    $\\phi(x) = [1, x, x^2, \\dots, x^{d-1}]$\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from numpy.polynomial.polynomial import polyvander\n",
    "def featurize_vandermonde(x, d, normalize = False):\n",
    "    A = polyvander(x, d-1)\n",
    "    for d_ in range(A.shape[1]):\n",
    "        if normalize:\n",
    "            A[:,d_] *=  np.sqrt(2*d_+1)\n",
    "    return A\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_poly_features(d):\n",
    "    n = 128\n",
    "    d_max = 20\n",
    "    x_type = 'uniform_random'\n",
    "    np.random.seed(7)\n",
    "    x_true = generate_x(x_type = 'grid', n=1000)\n",
    "    x_train = generate_x(x_type=x_type, n=n)\n",
    "    phi_train = featurize_vandermonde(x_train, d_max)\n",
    "    phi_true = featurize_vandermonde(x_true, d_max)\n",
    "\n",
    "    plt.plot(x_true, phi_true[:,d], linewidth = 0.5)\n",
    "    plt.scatter(x_train, phi_train[:,d], marker='o')\n",
    "    plt.ylim([-1.2,1.2])\n",
    "    plt.xlabel('x')\n",
    "    plt.ylabel('$\\phi(x)$')\n",
    "    plt.show()\n",
    "\n",
    "\n",
    "slider = widgets.IntSlider(\n",
    "    value=0,\n",
    "    min=0,\n",
    "    max=10,\n",
    "    step=1,\n",
    "    description='Feature # k:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "\n",
    "\n",
    "interactive_plot = interactive(plot_poly_features, d=slider)\n",
    "output = interactive_plot.children[-1]\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Fourier features\n",
    "We consider the d-dimensional real Fourier features given by:<br>\n",
    "    $\\phi(x) = [1, \\sin(\\pi x), \\cos(\\pi x), \\sin(2 \\pi x), \\cos(2\\pi x), \\dots,  \\sin (r \\pi x), \\cos(r \\pi x)]$, <br>\n",
    "    where $r = \\frac{d-1}{2}$.\n",
    "\n",
    "Note that by this convention we require $d$ to be an odd integer.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from numpy.polynomial.polynomial import polyvander\n",
    "def featurize_fourier(x, d, normalize = False):\n",
    "    assert (d-1) % 2 == 0, \"d must be odd\"\n",
    "    max_r = int((d-1)/2)\n",
    "    n = len(x)\n",
    "    A = np.zeros((n, d))\n",
    "    A[:,0] = 1\n",
    "    for d_ in range(1,max_r+1):\n",
    "        A[:,2*(d_-1)+1] =  np.sin(d_*x*np.pi)\n",
    "        A[:,2*(d_-1)+2] =  np.cos(d_*x*np.pi)\n",
    "\n",
    "    if normalize:\n",
    "        A[:,0] *= (1/np.sqrt(2))\n",
    "        A *= np.sqrt(2)\n",
    "    return A\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_fourier_features(x_type,d):\n",
    "    n = 128\n",
    "    d_max = 21\n",
    "    np.random.seed(7)\n",
    "    x_true = generate_x(x_type = 'grid', n=1000)\n",
    "    x_train = generate_x(x_type=x_type, n=n)\n",
    "    phi_train = featurize_fourier(x_train, d_max)\n",
    "    phi_true = featurize_fourier(x_true, d_max)\n",
    "\n",
    "    plt.plot(x_true, phi_true[:,d], linewidth = 0.5)\n",
    "    plt.scatter(x_train, phi_train[:,d], marker='o')\n",
    "    plt.ylim([-1.2,1.2])\n",
    "    plt.xlabel('x')\n",
    "    plt.ylabel('$\\phi(x)$')\n",
    "    plt.show()\n",
    "\n",
    "\n",
    "slider1 = widgets.RadioButtons(\n",
    "    options=['uniform_random', 'grid'],\n",
    "    description='x_type:',\n",
    "    disabled=False\n",
    ")\n",
    "\n",
    "\n",
    "slider2 = widgets.IntSlider(\n",
    "    value=0,\n",
    "    min=0,\n",
    "    max=20,\n",
    "    step=1,\n",
    "    description='Feature # k:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "\n",
    "\n",
    "interactive_plot = interactive(plot_fourier_features, d=slider2, x_type=slider1)\n",
    "output = interactive_plot.children[-1]\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def featurize(x, d, phi_type, normalize = False):\n",
    "    function_map = {'polynomial':featurize_vandermonde, 'fourier':featurize_fourier}\n",
    "    return function_map[phi_type](x,d,normalize)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part c): Linear Regression to learn the 1d-function in feature space\n",
    "\n",
    "Fill in the code in the functions solve_ls and solve_ridge to populate a variable \"coeffs\" that contains the minimizer $w^*$. We want you to learn how to use  sklearn. To check for correctness you can compare your solutions using sklearn to that obtained using the explicit formulas/numpy's least squares function.\n",
    "\n",
    "Then run cells in part c).  Report the learned coeffs $w^*$ and the training loss value from \"CELL LS_LINE\".\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "#### Least squares\n",
    "To learn the function we will perform linear regression in the lifted feature space, i.e. we learn a set of coefficients $w \\in \\mathbb{R}^d$ to minimize the least-squares loss:\n",
    "\n",
    "$\\ell(w) = \\frac{1}{n} \\| y - \\phi w \\|_2^2$.\n",
    "\n",
    "Assume that all training samples $x_i$ are unique.\n",
    "Here $\\phi$ has dimensions $n \\times d$. Depending on how $n$ and $d$ are related we are in one of the following regimes:\n",
    "1. $n > d$. This is the underparameterized regime and we have more samples than free parameters and there exists a unique solution that minimizes our loss, <br>\n",
    "$w^* = (\\phi^\\top \\phi)^{-1} \\phi^\\top y$.\n",
    "2. $n=d$. In this case the $\\phi$ matrix is square and invertible and the solution is given by,\n",
    "<br>\n",
    "$w^* = (\\phi)^{-1}  y$.\n",
    "3. $n < d$. This is the overparameterized regime and there are infinitely many solutions that acheive zero training loss. In this problem we consider the minimum $\\ell_2$ norm solution given by,<br>\n",
    "$w^* =  \\phi^\\top  ( \\phi \\phi^\\top)^{-1} y$.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#LABEL: CELL LR\n",
    "\n",
    "from sklearn.linear_model import LinearRegression\n",
    "def solve_ls(phi, y):\n",
    "\n",
    "    LR = LinearRegression(fit_intercept=False, normalize=False)\n",
    "\n",
    "    ### start c1 ###\n",
    "\n",
    "    ### end c1 ###\n",
    "\n",
    "    loss = np.mean((y- phi@coeffs)**2)\n",
    "    return coeffs, loss\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Ridge\n",
    "\n",
    "We can  add a regularizing penalty term to the least squares objective to perform ridge regression where we minimize the loss,\n",
    "\n",
    "$\\ell(w) = \\frac{1}{n} (\\| y - \\phi w \\|_2^2+ \\lambda \\| w \\|_2^2)$.\n",
    "\n",
    "In this case, irrespective of n and d, the minimizer is unique and is given by,\n",
    "\n",
    "$w^* =  \\phi^\\top  ( \\phi \\phi^\\top + \\lambda I_n)^{-1} y = (\\phi^\\top \\phi + \\lambda I_d)^{-1} \\phi^\\top y.$\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.linear_model import Ridge\n",
    "\n",
    "def solve_ridge(phi, y, lambda_ridge):\n",
    "\n",
    "    Rdg = Ridge(fit_intercept=False, normalize=False, alpha = lambda_ridge)\n",
    "\n",
    "    ### start c2 ###\n",
    "\n",
    "    ### end c2 ###\n",
    "\n",
    "    loss = np.mean((y- phi@coeffs)**2) + lambda_ridge*np.mean(coeffs**2)\n",
    "    return coeffs, loss\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Learning function using least squares\n",
    "Let us use least squares and try to learn the function $y = x$ from noisy training samples using polynomial features.\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#LABEL: CELL LS_LINE\n",
    "\n",
    "def get_params1():\n",
    "    np.random.seed(18)\n",
    "    n = 64\n",
    "    d = 2\n",
    "    awgn_std = 1e-1\n",
    "    x_type = 'uniform_random'\n",
    "    phi_type = 'polynomial'\n",
    "    f_type = 'x'\n",
    "\n",
    "    return n, d, awgn_std, x_type, phi_type, f_type\n",
    "\n",
    "n, d, awgn_std, x_type, phi_type, f_type = get_params1()\n",
    "\n",
    "x_train = generate_x(x_type=x_type, n=n)\n",
    "phi_train = featurize(x_train, d, phi_type)\n",
    "y_train = generate_y(x=x_train, f_type=f_type)\n",
    "y_train = add_awgn_noise(y_train, awgn_std)\n",
    "w, loss = solve_ls(phi_train , y_train)\n",
    "\n",
    "# print(\"w:\", w)\n",
    "# print(\"loss:\", loss)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part d)  Visualizing the learned function\n",
    "\n",
    "Fill in the code in function get_plot_data to populate the variable \"y_plot_pred\" that contains the function values predicted by our learned function. Then run the cells in part d) to visualize the learned function plotted alongside the true function.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualizing learned function vs true function\n",
    "\n",
    "Now that we have a set of coefficients that determines our learned function we would like to plot the learned function alongside the true function to visually determine how good our approximation is\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, d, awgn_std, x_type, phi_type, f_type = get_params1()\n",
    "\n",
    "def get_plot_data( f_type, phi_type, d, w, n_plot = 1000):\n",
    "    x_plot= generate_x(x_type = 'grid', n=n_plot)\n",
    "    y_plot_true = generate_y(x=x_plot, f_type=f_type)\n",
    "\n",
    "    #Computing predictions via learned function\n",
    "\n",
    "    phi_plot = featurize(x_plot, d, phi_type)\n",
    "\n",
    "    ### start d ###\n",
    "\n",
    "    ### end d ###\n",
    "\n",
    "    return x_plot, y_plot_true, y_plot_pred\n",
    "\n",
    "\n",
    "x_plot, y_plot_true, y_plot_pred = get_plot_data(f_type, phi_type, d, w, n_plot = 1000)\n",
    "plt.plot(x_plot, y_plot_pred, 'o-', ms=2, label = 'Learned function')\n",
    "plt.plot(x_plot, y_plot_true, label = 'True function')\n",
    "plt.scatter(x_train, y_train, marker='o', s=10, label = 'Training samples')\n",
    "plt.legend()\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part e)  Calculating goodness of approximation/fit\n",
    "\n",
    "Fill in code in function get_fit_mse to populate the variable \"fit_error\" that contains the fit mean squared-error an empirical approximation to the distance between the true function and the learned function. Run cells in part e) and report the fit mean squared error.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "To quantify the goodness of approximation/fit we would like to measure the \"distance\" between the learned function and the true function. Since we are interested in the functions in the range (-1,1) one way to measure distance is to compute the integral,\n",
    "$d(f_{true}, f_{learned}) = \\frac{1}{2} \\int_{-1}^{1} (f_{true} (x) - f_{learned}(x))^2 dx$\n",
    "\n",
    "We can approximate this integral by taking an empirical mean of the squared difference between the function values of the true and learned function over a batch. We construct this batch by sampling uniform randomly from the interval $[-1,1]$.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#LABEL: CELL FIT MSE\n",
    "n, d, awgn_std, x_type, phi_type, f_type = get_params1()\n",
    "\n",
    "def get_fit_mse( f_type, phi_type, d, w, n_fit = 1000):\n",
    "    x_fit= generate_x(x_type = 'uniform_random', n=n_fit)\n",
    "    y_fit_true = generate_y(x=x_fit, f_type=f_type)\n",
    "\n",
    "    #Computing predictions via learned function\n",
    "\n",
    "    phi_fit = featurize(x_fit, d, phi_type)\n",
    "    y_fit_pred = phi_fit @ w\n",
    "\n",
    "    ### start e ###\n",
    "\n",
    "    ### end e ###\n",
    "    return fit_error\n",
    "\n",
    "\n",
    "fit_mse = get_fit_mse(f_type, phi_type, d, w, n_fit = 1000)\n",
    "\n",
    "# print(\"Fit MSE: \", fit_mse)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part f): Choice of feature type matters\n",
    "Run the cells in part f) and explore the effect of number of features $d$, and the choice of feature type in learning the functions $y = x$, and $y = cos(2\\pi x)$. What do you observe?\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "While using both fourier and polynomial features we can approximate any function in an interval, the number of features that we need to consider varies based on the function we are trying to approximate. Some functions are easy to approximate using polynomial features while others are easier to approximate using Fourier features.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def run1(n, d, x_type, f_type, phi_type, seed = 1, awgn_std = 0, lambda_ridge = 0):\n",
    "\n",
    "\n",
    "    #Set randomness\n",
    "    np.random.seed(seed)\n",
    "\n",
    "    x_train = generate_x(x_type=x_type, n=n)\n",
    "    phi_train = featurize(x_train, d, phi_type)\n",
    "    y_train = generate_y(x=x_train, f_type=f_type)\n",
    "\n",
    "    if awgn_std != 0:\n",
    "        y_train = add_awgn_noise(y_train, awgn_std)\n",
    "\n",
    "    if lambda_ridge == 0:\n",
    "\n",
    "        #Write one line of code calling the function to solve least squares that returns w, loss\n",
    "        ### start f1 ###\n",
    "\n",
    "        ### end f1 ###\n",
    "    else:\n",
    "        #Write one line of code calling the function to solve ridge regression that returns w, loss\n",
    "\n",
    "        ### start f2 ###\n",
    "\n",
    "        ### end f2 ###\n",
    "\n",
    "    return x_train, y_train, w, loss\n",
    "\n",
    "\n",
    "def visualize1(x_train, y_train, f_type, phi_type, d, w, loss, n_plot = 1000, n_fit = 1000):\n",
    "    x_plot, y_plot_true, y_plot_pred = get_plot_data(f_type, phi_type, d, w, n_plot)\n",
    "    plt.plot(x_plot, y_plot_true, label = 'True function')\n",
    "    plt.scatter(x_train, y_train, marker='o', s=20, label = 'Training samples')\n",
    "    plt.plot(x_plot, y_plot_pred, 'o-', ms=2, label = 'Learned function')\n",
    "\n",
    "\n",
    "    fit_mse = get_fit_mse(f_type, phi_type, d, w, n_fit = n_fit)\n",
    "    plt.title(\"Train loss:\" + str(\"{:.2e}\".format(loss)) + \", Fit MSE: \" + str(\"{:.2e}\".format(fit_mse)))\n",
    "    plt.ylim([-1.5, 1.5])\n",
    "    plt.xlabel('x')\n",
    "    plt.ylabel('y')\n",
    "    plt.legend(bbox_to_anchor = (1.03, 0.97))\n",
    "    plt.show()\n",
    "\n",
    "#     plt.plot(w, 'o')\n",
    "    markerlines, stemlines,  baseline = plt.stem(np.arange(d), w, 'b', 'o',  use_line_collection=True)\n",
    "    plt.setp(stemlines, 'color', plt.getp(markerlines,'color'))\n",
    "    plt.xlabel('feature #(k)')\n",
    "    plt.ylabel('weight')\n",
    "    plt.show()\n",
    "\n",
    "\n",
    "    return fit_mse\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_feature_type_matters(n, n_plot, n_fit, x_type, f_type, phi_type, seed, awgn_std, lambda_ridge, d):\n",
    "    x_train, y_train, w, loss = run1(n, d, x_type, f_type, phi_type, awgn_std = awgn_std, lambda_ridge = lambda_ridge, seed = seed)\n",
    "    fit_mse = visualize1(x_train, y_train, f_type, phi_type, d, w, loss, n_plot , n_fit)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#CELL FTM 1\n",
    "def get_params2():\n",
    "    n = 64\n",
    "    n_plot = 1000\n",
    "    n_fit = 10000\n",
    "    x_type = 'uniform_random'\n",
    "    seed = 1\n",
    "    awgn_std = 0\n",
    "    lambda_ridge = 0\n",
    "    return n, n_plot, n_fit, x_type, seed, awgn_std, lambda_ridge\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, n_plot, n_fit, x_type, seed, awgn_std, lambda_ridge = get_params2()\n",
    "slider1 = widgets.IntSlider(\n",
    "    value=1,\n",
    "    min=1,\n",
    "    max=25,\n",
    "    step=1,\n",
    "    description='d:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "\n",
    "f_type = 'x'\n",
    "phi_type = 'polynomial'\n",
    "print(\"Polynomial features\")\n",
    "\n",
    "interactive_plot =interactive(plot_feature_type_matters,n = fixed(n), n_plot = fixed(n_plot), n_fit = fixed(n_fit),\n",
    "                              x_type = fixed(x_type),f_type = fixed(f_type), phi_type = fixed(phi_type),seed = fixed(seed),\n",
    "                              awgn_std= fixed(awgn_std), lambda_ridge = fixed(lambda_ridge), d = slider1)\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, n_plot, n_fit, x_type, seed, awgn_std, lambda_ridge = get_params2()\n",
    "\n",
    "slider2 = widgets.IntSlider(\n",
    "    value=1,\n",
    "    min=1,\n",
    "    max=25,\n",
    "    step=2,\n",
    "    description='d:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "\n",
    "f_type = 'x'\n",
    "phi_type = 'fourier'\n",
    "print(\"Fourier features\")\n",
    "\n",
    "interactive_plot =interactive(plot_feature_type_matters,n = fixed(n), n_plot = fixed(n_plot), n_fit = fixed(n_fit),\n",
    "                              x_type = fixed(x_type),f_type = fixed(f_type), phi_type = fixed(phi_type),seed = fixed(seed),\n",
    "                              awgn_std= fixed(awgn_std), lambda_ridge = fixed(lambda_ridge), d = slider2)\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, n_plot, n_fit, x_type, seed, awgn_std, lambda_ridge = get_params2()\n",
    "\n",
    "slider1 = widgets.IntSlider(\n",
    "    value=1,\n",
    "    min=1,\n",
    "    max=25,\n",
    "    step=1,\n",
    "    description='d:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "\n",
    "f_type = 'cos2'\n",
    "phi_type = 'polynomial'\n",
    "print(\"Polynomial features\")\n",
    "interactive_plot =interactive(plot_feature_type_matters,n = fixed(n), n_plot = fixed(n_plot), n_fit = fixed(n_fit),\n",
    "                              x_type = fixed(x_type),f_type = fixed(f_type), phi_type = fixed(phi_type),seed = fixed(seed),\n",
    "                              awgn_std= fixed(awgn_std), lambda_ridge = fixed(lambda_ridge), d = slider1)\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, n_plot, n_fit, x_type, seed, awgn_std, lambda_ridge = get_params2()\n",
    "\n",
    "slider2 = widgets.IntSlider(\n",
    "    value=1,\n",
    "    min=1,\n",
    "    max=25,\n",
    "    step=2,\n",
    "    description='d:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "\n",
    "f_type = 'cos2'\n",
    "phi_type = 'fourier'\n",
    "print(\"Fourier features\")\n",
    "\n",
    "interactive_plot =interactive(plot_feature_type_matters,n = fixed(n), n_plot = fixed(n_plot), n_fit = fixed(n_fit),\n",
    "                              x_type = fixed(x_type),f_type = fixed(f_type), phi_type = fixed(phi_type),seed = fixed(seed),\n",
    "                              awgn_std= fixed(awgn_std), lambda_ridge = fixed(lambda_ridge), d = slider2)\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Measuring the approxmiation error vs d\n",
    "def plot_errors(n, d_range, n_fit, x_type, seed, awgn_std, lambda_ridge, phi_type, f_type, title = None, ylim = -1):\n",
    "    train_loss_array = []\n",
    "    fit_mse_array = []\n",
    "    for d in d_range:\n",
    "        x_train, y_train, w, loss = run1(n, d, x_type, f_type, phi_type, awgn_std=awgn_std, lambda_ridge = lambda_ridge, seed=seed)\n",
    "        fit_mse = get_fit_mse(f_type, phi_type, d, w, n_fit = n_fit)\n",
    "        train_loss_array.append(loss)\n",
    "        fit_mse_array.append(fit_mse)\n",
    "\n",
    "    plt.plot(d_range, train_loss_array, 'o-', label = 'Train loss')\n",
    "    plt.plot(d_range, fit_mse_array, 'o-', label = 'Fit MSE')\n",
    "    plt.yscale('log')\n",
    "    plt.ylabel('Loss/Error')\n",
    "    plt.xlabel('d')\n",
    "    if title is not None:\n",
    "        plt.title(title)\n",
    "    plt.legend(bbox_to_anchor = (1.33, 0.97))\n",
    "\n",
    "\n",
    "    if ylim != -1:\n",
    "        plt.ylim(ylim)\n",
    "    plt.show()\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 64\n",
    "n_fit = 10000\n",
    "x_type = 'uniform_random'\n",
    "seed = 1\n",
    "awgn_std = 0\n",
    "lambda_ridge = 0\n",
    "\n",
    "phi_type = 'polynomial'\n",
    "for f_type in ['x', 'cos2']:\n",
    "    d_range = np.arange(1,25,1)\n",
    "    plot_errors(n, d_range, n_fit, x_type, seed, awgn_std, lambda_ridge, phi_type, f_type, title = 'f_type:' + str(f_type) + ', phi_type: '  + str(phi_type))\n",
    "\n",
    "\n",
    "phi_type = 'fourier'\n",
    "for f_type in ['x', 'cos2']:\n",
    "    d_range = np.arange(1,25,2)\n",
    "    plot_errors(n, d_range, n_fit, x_type, seed, awgn_std, lambda_ridge, phi_type, f_type, title = 'f_type:' + str(f_type) + ', phi_type: '  + str(phi_type))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part g): Effect of noise\n",
    "Run the cells in part g) and explore the effect of noise on our learning process. What do you observe as d grows to n? Does one particualr combination of input sampling/ feature type do better than the others?\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "Next we will see the effect of noise in our learning method. For this purpose we will consider the three kinds of input data/feature combinations:\n",
    "1. uniform randomly sampled x, polynomial features\n",
    "2. uniform randomly sampled x, fourier features\n",
    "3. evenly spaced x, fourier features\n",
    "\n",
    "Further in this case, we will see what happens when d grows all the way to n.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#CELL EON 1\n",
    "def get_params3():\n",
    "    n = 65\n",
    "    n_plot = 1000\n",
    "    n_fit = 10000\n",
    "    seed = 1\n",
    "    awgn_std = 1e-1\n",
    "    lambda_ridge = 0\n",
    "    return n, n_plot, n_fit, seed, awgn_std, lambda_ridge\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, n_plot, n_fit, seed, awgn_std, lambda_ridge = get_params3()\n",
    "slider1 = widgets.IntSlider(\n",
    "    value=1,\n",
    "    min=1,\n",
    "    max=65,\n",
    "    step=1,\n",
    "    description='d:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "x_type = 'uniform_random'\n",
    "\n",
    "f_type = 'cos2'\n",
    "phi_type = 'polynomial'\n",
    "print(\"Uniform random x, Polynomial features\")\n",
    "interactive_plot =interactive(plot_feature_type_matters,n = fixed(n), n_plot = fixed(n_plot), n_fit = fixed(n_fit),\n",
    "                              x_type = fixed(x_type),f_type = fixed(f_type), phi_type = fixed(phi_type),seed = fixed(seed),\n",
    "                              awgn_std= fixed(awgn_std), lambda_ridge = fixed(lambda_ridge), d = slider1)\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, n_plot, n_fit, seed, awgn_std, lambda_ridge = get_params3()\n",
    "slider2 = widgets.IntSlider(\n",
    "    value=1,\n",
    "    min=1,\n",
    "    max=65,\n",
    "    step=2,\n",
    "    description='d:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "x_type = 'uniform_random'\n",
    "\n",
    "f_type = 'cos2'\n",
    "phi_type = 'fourier'\n",
    "print(\"Uniform random x, Fourier features\")\n",
    "\n",
    "interactive_plot =interactive(plot_feature_type_matters,n = fixed(n), n_plot = fixed(n_plot), n_fit = fixed(n_fit),\n",
    "                              x_type = fixed(x_type),f_type = fixed(f_type), phi_type = fixed(phi_type),seed = fixed(seed),\n",
    "                              awgn_std= fixed(awgn_std), lambda_ridge = fixed(lambda_ridge), d = slider2)\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, n_plot, n_fit, seed, awgn_std, lambda_ridge = get_params3()\n",
    "slider2 = widgets.IntSlider(\n",
    "    value=1,\n",
    "    min=1,\n",
    "    max=65,\n",
    "    step=2,\n",
    "    description='d:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "x_type = 'grid'\n",
    "\n",
    "f_type = 'cos2'\n",
    "phi_type = 'fourier'\n",
    "print(\"Evenly spaced x, Fourier features\")\n",
    "\n",
    "interactive_plot =interactive(plot_feature_type_matters,n = fixed(n), n_plot = fixed(n_plot), n_fit = fixed(n_fit),\n",
    "                              x_type = fixed(x_type),f_type = fixed(f_type), phi_type = fixed(phi_type),seed = fixed(seed),\n",
    "                              awgn_std= fixed(awgn_std), lambda_ridge = fixed(lambda_ridge), d = slider2)\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 65\n",
    "n_fit = 10000\n",
    "\n",
    "seed = 1\n",
    "awgn_std = 1e-1\n",
    "lambda_ridge = 0\n",
    "f_type = 'cos2'\n",
    "\n",
    "x_type = 'uniform_random'\n",
    "phi_type = 'polynomial'\n",
    "d_range = np.arange(1,65,1)\n",
    "plot_errors(n, d_range, n_fit, x_type, seed, awgn_std, lambda_ridge, phi_type, f_type, title = 'x_type: ' + str(x_type) + ', phi_type: ' + str(phi_type))\n",
    "\n",
    "x_type = 'uniform_random'\n",
    "phi_type = 'fourier'\n",
    "d_range = np.arange(1,65,2)\n",
    "plot_errors(n, d_range, n_fit, x_type, seed, awgn_std, lambda_ridge, phi_type, f_type, title = 'x_type: ' + str(x_type) + ', phi_type: ' + str(phi_type))\n",
    "\n",
    "x_type = 'grid'\n",
    "phi_type = 'fourier'\n",
    "d_range = np.arange(1,65,2)\n",
    "plot_errors(n, d_range, n_fit, x_type, seed, awgn_std, lambda_ridge, phi_type, f_type, title = 'x_type: ' + str(x_type) + ', phi_type: ' + str(phi_type))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part h) conditioning of $\\phi^T \\phi$\n",
    "Run the cells in part h) and visualize the eigenvalue spectrum of $\\phi^T \\phi$ we see that the matrix can get highly ill conditioned depending on the input sampling and feature type.  Comment on the shape of the curve for poly random. Do you think thats really what the curve should look like? If not can you think of a reason why this is happening?\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_eig_values(n, d, seed, lambda_ridge, shadow = False):\n",
    "    np.random.seed(seed)\n",
    "    x_type_phi_type_pairs = [('uniform_random', 'polynomial'), ('uniform_random', 'fourier'), ('grid', 'fourier')]\n",
    "\n",
    "    colors = ['blue', 'green', 'orange']\n",
    "\n",
    "    for k, (x_type, phi_type) in enumerate(x_type_phi_type_pairs):\n",
    "        x_train = generate_x(x_type=x_type, n=n)\n",
    "        phi_train = featurize(x_train, d, phi_type)\n",
    "        eig_vals,_ = np.linalg.eig(phi_train.T @ phi_train + lambda_ridge*np.eye(d))\n",
    "\n",
    "#         print(eig_vals)\n",
    "\n",
    "        eig_vals = np.sort(np.abs(eig_vals))[::-1]\n",
    "        plt.plot(eig_vals, 'o-', c = colors[k], label = 'x_type: ' + str(x_type) + ', phi_type: ' + str(phi_type))\n",
    "        if shadow is True and lambda_ridge != 0:\n",
    "\n",
    "            eig_vals_shadow,_ = np.linalg.eig(phi_train.T @ phi_train )\n",
    "            eig_vals_shadow = np.sort(np.abs(eig_vals_shadow))[::-1]\n",
    "\n",
    "            plt.plot(eig_vals_shadow, '-', c = colors[k], label = '$\\lambda=0$, x_type: ' + str(x_type) + ', phi_type: ' + str(phi_type), alpha = 0.5)\n",
    "\n",
    "        #         print(eig_vals)\n",
    "\n",
    "\n",
    "    plt.legend(bbox_to_anchor = (1.73, 0.97))\n",
    "    plt.yscale('log')\n",
    "    plt.ylim(1e-20, 1e4)\n",
    "    plt.xlim([-1, n+1])\n",
    "\n",
    "    plt.show()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 65\n",
    "seed = 1\n",
    "awgn_std = 1e-1\n",
    "lambda_ridge = 0\n",
    "\n",
    "\n",
    "slider = widgets.IntSlider(\n",
    "    value=11,\n",
    "    min=1,\n",
    "    max=65,\n",
    "    step=2,\n",
    "    description='d:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "\n",
    "interactive_plot =interactive(plot_eig_values,n = fixed(n), d = slider, seed = fixed(seed),\n",
    "                              lambda_ridge = fixed(lambda_ridge), shadow =fixed(False))\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part i) Ridge regression to combat noise\n",
    "Run the cells in part i) and explore the effect of using ridge regression instead of least squares. Can we use ridge regression to combat the effect of noise and allow us to get good approximation even using $d$ close to $n$? Play around. Can you adjust $\\lambda$ to combat the effect of noise at d = n?\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#CELL EOR 1\n",
    "def get_params4():\n",
    "    n = 65\n",
    "    n_plot = 1000\n",
    "    n_fit = 10000\n",
    "    seed = 1\n",
    "    awgn_std = 1e-1\n",
    "    return n, n_plot, n_fit, seed, awgn_std\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, n_plot, n_fit, seed, awgn_std = get_params4()\n",
    "\n",
    "slider1 = widgets.IntSlider(\n",
    "    value=65,\n",
    "    min=1,\n",
    "    max=65,\n",
    "    step=1,\n",
    "    description='d:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "\n",
    "slider3= widgets.FloatLogSlider(\n",
    "    value=-25,\n",
    "    base=10,\n",
    "    min=-25, # max exponent of base\n",
    "    max=2, # min exponent of base\n",
    "    step=1, # exponent step\n",
    "    description='$\\lambda$',\n",
    "    continuous_update= False\n",
    ")\n",
    "\n",
    "x_type = 'uniform_random'\n",
    "\n",
    "f_type = 'cos2'\n",
    "phi_type = 'polynomial'\n",
    "print(\"Uniform random x, Polynomial features\")\n",
    "interactive_plot =interactive(plot_feature_type_matters,n = fixed(n), n_plot = fixed(n_plot), n_fit = fixed(n_fit),\n",
    "                              x_type = fixed(x_type),f_type = fixed(f_type), phi_type = fixed(phi_type),seed = fixed(seed),\n",
    "                              awgn_std= fixed(awgn_std), lambda_ridge = slider3, d = slider1)\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, n_plot, n_fit, seed, awgn_std = get_params4()\n",
    "slider2 = widgets.IntSlider(\n",
    "    value=65,\n",
    "    min=1,\n",
    "    max=65,\n",
    "    step=2,\n",
    "    description='d:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "\n",
    "slider3= widgets.FloatLogSlider(\n",
    "    value=-25,\n",
    "    base=10,\n",
    "    min=-25, # max exponent of base\n",
    "    max=2, # min exponent of base\n",
    "    step=1, # exponent step\n",
    "    description='$\\lambda$',\n",
    "    continuous_update= False\n",
    ")\n",
    "\n",
    "x_type = 'uniform_random'\n",
    "\n",
    "f_type = 'cos2'\n",
    "phi_type = 'fourier'\n",
    "print(\"Uniform random x, Fourier features\")\n",
    "\n",
    "interactive_plot =interactive(plot_feature_type_matters,n = fixed(n), n_plot = fixed(n_plot), n_fit = fixed(n_fit),\n",
    "                              x_type = fixed(x_type),f_type = fixed(f_type), phi_type = fixed(phi_type),seed = fixed(seed),\n",
    "                              awgn_std= fixed(awgn_std), lambda_ridge = slider3, d = slider2)\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, n_plot, n_fit, seed, awgn_std = get_params4()\n",
    "slider2 = widgets.IntSlider(\n",
    "    value=65,\n",
    "    min=1,\n",
    "    max=65,\n",
    "    step=2,\n",
    "    description='d:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "\n",
    "slider3= widgets.FloatLogSlider(\n",
    "    value=-25,\n",
    "    base=10,\n",
    "    min=-25, # max exponent of base\n",
    "    max=2, # min exponent of base\n",
    "    step=1, # exponent step\n",
    "    description='$\\lambda$',\n",
    "    continuous_update= False\n",
    ")\n",
    "\n",
    "x_type = 'grid'\n",
    "\n",
    "f_type = 'cos2'\n",
    "phi_type = 'fourier'\n",
    "print(\"Evenly spaced x, Fourier features\")\n",
    "\n",
    "interactive_plot =interactive(plot_feature_type_matters,n = fixed(n), n_plot = fixed(n_plot), n_fit = fixed(n_fit),\n",
    "                              x_type = fixed(x_type),f_type = fixed(f_type), phi_type = fixed(phi_type),seed = fixed(seed),\n",
    "                              awgn_std= fixed(awgn_std), lambda_ridge = slider3, d = slider2)\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Plotting the fit errors\n",
    "def get_params5():\n",
    "    n = 65\n",
    "    n_fit = 10000\n",
    "\n",
    "    seed = 1\n",
    "    awgn_std = 1e-1\n",
    "    # lambda_ridge = 1e-1\n",
    "    f_type = 'cos2'\n",
    "    return n, n_fit, seed, awgn_std, f_type\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, n_fit, seed, awgn_std, f_type = get_params5()\n",
    "x_type = 'uniform_random'\n",
    "phi_type = 'polynomial'\n",
    "slider3= widgets.FloatLogSlider(\n",
    "    value=1e-1,\n",
    "    base=10,\n",
    "    min=-6, # max exponent of base\n",
    "    max=2, # min exponent of base\n",
    "    step=1, # exponent step\n",
    "    description='$\\lambda$',\n",
    "    continuous_update= False\n",
    ")\n",
    "d_range = np.arange(1,65,1)\n",
    "\n",
    "interactive_plot =interactive(plot_errors, ylim = fixed([1e-4, 1e5]), n = fixed(n), n_fit = fixed(n_fit),  d_range = fixed(d_range),\n",
    "                              x_type = fixed(x_type),f_type = fixed(f_type), phi_type = fixed(phi_type),seed = fixed(seed),\n",
    "                              awgn_std= fixed(awgn_std), title = fixed('x_type: ' + str(x_type) + ', phi_type: ' + str(phi_type)),  lambda_ridge = slider3, d = slider2)\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, n_fit, seed, awgn_std, f_type = get_params5()\n",
    "\n",
    "x_type = 'uniform_random'\n",
    "phi_type = 'fourier'\n",
    "slider3= widgets.FloatLogSlider(\n",
    "    value=1e-1,\n",
    "    base=10,\n",
    "    min=-6, # max exponent of base\n",
    "    max=2, # min exponent of base\n",
    "    step=1, # exponent step\n",
    "    description='$\\lambda$',\n",
    "    continuous_update= False\n",
    ")\n",
    "d_range = np.arange(1,65,2)\n",
    "\n",
    "interactive_plot =interactive(plot_errors, ylim = fixed([1e-4, 1e5]), n = fixed(n), n_fit = fixed(n_fit),  d_range = fixed(d_range),\n",
    "                              x_type = fixed(x_type),f_type = fixed(f_type), phi_type = fixed(phi_type),seed = fixed(seed),\n",
    "                              awgn_std= fixed(awgn_std), title = fixed('x_type: ' + str(x_type) + ', phi_type: ' + str(phi_type)),  lambda_ridge = slider3, d = slider2)\n",
    "interactive_plot\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n, n_fit, seed, awgn_std, f_type = get_params5()\n",
    "\n",
    "x_type = 'grid'\n",
    "phi_type = 'fourier'\n",
    "slider3= widgets.FloatLogSlider(\n",
    "    value=1e-1,\n",
    "    base=10,\n",
    "    min=-6, # max exponent of base\n",
    "    max=2, # min exponent of base\n",
    "    step=1, # exponent step\n",
    "    description='$\\lambda$',\n",
    "    continuous_update= False\n",
    ")\n",
    "d_range = np.arange(1,65,2)\n",
    "\n",
    "interactive_plot =interactive(plot_errors, ylim = fixed([1e-4, 1e5]), n = fixed(n), n_fit = fixed(n_fit),  d_range = fixed(d_range),\n",
    "                              x_type = fixed(x_type),f_type = fixed(f_type), phi_type = fixed(phi_type),seed = fixed(seed),\n",
    "                              awgn_std= fixed(awgn_std), title = fixed('x_type: ' + str(x_type) + ', phi_type: ' + str(phi_type)),  lambda_ridge = slider3, d = slider2)\n",
    "interactive_plot\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part j) conditioning of $\\phi^T \\phi  + \\lambda I$\n",
    "Run the cells in part j) and observe the eigenvalue spectrum. What is the effect of $\\lambda$ on the shape of the eigenvalue curve?\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 65\n",
    "d = 65\n",
    "seed = 1\n",
    "awgn_std = 1e-1\n",
    "lambda_ridge = 0\n",
    "slider3= widgets.FloatLogSlider(\n",
    "    value=1e-6,\n",
    "    base=10,\n",
    "    min=-20, # max exponent of base\n",
    "    max=2, # min exponent of base\n",
    "    step=1, # exponent step\n",
    "    description='$\\lambda$',\n",
    "    continuous_update= False\n",
    ")\n",
    "interactive_plot =interactive(plot_eig_values,n = fixed(n), d = fixed(d), seed = fixed(seed), lambda_ridge = slider3, shadow = fixed(True))\n",
    "interactive_plot\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
